#All Figures for Brother and Burden 

rm(list = ls())

#Clear the console
cat("\014")

library(foreign)
library(car)
library(gplots)
library(ggplot2)
library(plyr)
library(grid)

mydata <- read.dta("final data.dta")
attach(mydata)

### define MULTIPLOT function

multiplot <- function(..., plotlist=NULL, file, cols=1, layout=NULL) {
  require(grid)

  # Make a list from the ... arguments and plotlist
  plots <- c(list(...), plotlist)

  numPlots = length(plots)

  # If layout is NULL, then use 'cols' to determine layout
  if (is.null(layout)) {
    # Make the panel
    # ncol: Number of columns of plots
    # nrow: Number of rows needed, calculated from # of cols
    layout <- matrix(seq(1, cols * ceiling(numPlots/cols)),
                    ncol = cols, nrow = ceiling(numPlots/cols))
  }

 if (numPlots==1) {
    print(plots[[1]])

  } else {
    # Set up the page
    grid.newpage()
    pushViewport(viewport(layout = grid.layout(nrow(layout), ncol(layout))))

    # Make each plot, in the correct location
    for (i in 1:numPlots) {
      # Get the i,j matrix positions of the regions that contain this subplot
      matchidx <- as.data.frame(which(layout == i, arr.ind = TRUE))

      print(plots[[i]], vp = viewport(layout.pos.row = matchidx$row,
                                      layout.pos.col = matchidx$col))
    }
  }
}




####################### FIGURE 1: MAIN RESULTS


### Figure 1.1

control.stay <- shouldstay[treatment==1]
treatment_Sunni.stay <- shouldstay[treatment==2]
treatment_Muslim.stay <- shouldstay[treatment==3]
treatment_Cost.stay <- shouldstay[treatment==4]
treatment_Sunni_Cost.stay <- shouldstay[treatment==5]
treatment_Muslim_Cost.stay <- shouldstay[treatment==6]



groups <- seq(1,6,1)

groups.stay <- cbind(control.stay, treatment_Sunni.stay, treatment_Muslim.stay, treatment_Cost.stay, treatment_Sunni_Cost.stay, treatment_Muslim_Cost.stay)

means1 <- colMeans(groups.stay, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}



ses1 <- apply(groups.stay, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1


condition <- as.factor(c(1:6))
levels(condition) <- c("Control", "Sunni", "Muslim", "Cost", 
                       "Sunni/Cost", "Muslim/Cost")

data1 <- as.data.frame(cbind(condition, means1, ses1, uis1, lis1),
                       use.value.labels= T)

P1 <- ggplot(data1, aes(x=condition, y=means1)) +      
  geom_errorbar(aes(ymin=lis1, ymax=uis1),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(shape=21, size=3, fill = "black")+
  ylab("Mean response") +
  ggtitle("Figure 1.1. Refugees Can Stay \n") +
  coord_cartesian(ylim=c(0, 1))+ 
  scale_x_discrete(limits = c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))+
  theme(axis.title.x = element_blank()) 




### Figure 1.2

control.spend <- spend_ref[treatment==1]
treatment_Sunni.spend <- spend_ref[treatment==2]
treatment_Muslim.spend <- spend_ref[treatment==3]
treatment_Cost.spend <- spend_ref[treatment==4]
treatment_Sunni_Cost.spend <- spend_ref[treatment==5]
treatment_Muslim_Cost.spend <- spend_ref[treatment==6]


groups2 <- seq(1,6,1)

groups.spend <- cbind(control.spend, treatment_Sunni.spend, treatment_Muslim.spend, treatment_Cost.spend, treatment_Sunni_Cost.spend, treatment_Muslim_Cost.spend)

means2 <- colMeans(groups.spend, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses2 <- apply(groups.spend, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


data2 <- as.data.frame(cbind(condition, means2, ses2, uis2, lis2))


P2 <- ggplot(data2, aes(x=condition, y=means2)) +      
  geom_errorbar(aes(ymin=lis2, ymax=uis2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(shape=21, size=3, fill = "black")+
  ylab("Mean response") +
  ggtitle("Figure 1.2. Support for \nGovernment Spending") +
  coord_cartesian(ylim=c(1, 3))+  
  scale_x_discrete(limits = c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))+
  theme(axis.title.x = element_blank()) 




### Figure 1.3

control.trust <- trust_ref[treatment==1]
treatment_Sunni.trust <- trust_ref[treatment==2]
treatment_Muslim.trust <- trust_ref[treatment==3]
treatment_Cost.trust <- trust_ref[treatment==4]
treatment_Sunni_Cost.trust <- trust_ref[treatment==5]
treatment_Muslim_Cost.trust <- trust_ref[treatment==6]


groups3 <- seq(1,6,1)

groups.trust <- cbind(control.trust, treatment_Sunni.trust, treatment_Muslim.trust, treatment_Cost.trust, treatment_Sunni_Cost.trust, treatment_Muslim_Cost.trust)

means3 <- colMeans(groups.trust, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3 <- apply(groups.trust, 2, se.mean)

uis3 <- means3 + 1.96*ses3
lis3 <- means3 - 1.96*ses3

data3 <- as.data.frame(cbind(condition, means3, ses3, uis3, lis3))


P3 <- ggplot(data3, aes(x=condition, y=means3)) +      
  geom_errorbar(aes(ymin=lis3, ymax=uis3),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(shape=21, size=3, fill = "black")+
  ylab("Mean response") +
  ggtitle("Figure 1.3. Trust") +
  coord_cartesian(ylim=c(1, 3))+  
  scale_x_discrete(limits = c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))+
  theme(axis.title.x = element_blank()) 



### Figure 1.4


control.n <- neighborhood_ref[treatment==1]
treatment_Sunni.n <- neighborhood_ref[treatment==2]
treatment_Muslim.n <- neighborhood_ref[treatment==3]
treatment_Cost.n <- neighborhood_ref[treatment==4]
treatment_Sunni_Cost.n <- neighborhood_ref[treatment==5]
treatment_Muslim_Cost.n <- neighborhood_ref[treatment==6]



groups4 <- seq(1,6,1)

groups.n <- cbind(control.n, treatment_Sunni.n, treatment_Muslim.n, treatment_Cost.n, treatment_Sunni_Cost.n, treatment_Muslim_Cost.n)

means4 <- colMeans(groups.n, na.rm=TRUE)

se.mean4 <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses4 <- apply(groups.n, 2, se.mean4)

uis4 <- means4 + 1.96*ses4
lis4 <- means4 - 1.96*ses4

data4 <- as.data.frame(cbind(condition, means4, ses4, uis4, lis4))


P4 <- ggplot(data4, aes(x=condition, y=means4)) +      
  geom_errorbar(aes(ymin=lis4, ymax=uis4),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(shape=21, size=3, fill = "black")+
  ylab("Mean response") +
  ggtitle("Figure 1.4. Acceptance of the Refugees") +
  coord_cartesian(ylim=c(1, 3))+  
  scale_x_discrete(limits = c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))+
  theme(axis.title.x = element_blank()) 



### Figure 1.5

control.donation <- donation[treatment==1]
treatment_Sunni.donation <- donation[treatment==2]
treatment_Muslim.donation <- donation[treatment==3]
treatment_Cost.donation <- donation[treatment==4]
treatment_Sunni_Cost.donation <- donation[treatment==5]
treatment_Muslim_Cost.donation <- donation[treatment==6]


groups5 <- seq(1,6,1)

groups.donation <- cbind(control.donation, treatment_Sunni.donation, treatment_Muslim.donation, treatment_Cost.donation, treatment_Sunni_Cost.donation, treatment_Muslim_Cost.donation)

means5 <- colMeans(groups.donation, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses5 <- apply(groups.donation, 2, se.mean)

uis5 <- means5 + 1.96*ses5
lis5 <- means5 - 1.96*ses5


data5 <- as.data.frame(cbind(condition, means5, ses5, uis5, lis5))


P5<- ggplot(data5, aes(x=condition, y=means5)) +      
  geom_errorbar(aes(ymin=lis5, ymax=uis5),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(shape=21, size=3, fill = "black")+
  ylab("Mean donation") +
  ggtitle("Figure 1.5. Average Donation") +
  coord_cartesian(ylim=c(0, 4))+  
  scale_x_discrete(limits = c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))+
  theme(axis.title.x = element_blank()) 




### Figure 1.6

control.index <- support_general[treatment==1]
treatment_Sunni.index <- support_general[treatment==2]
treatment_Muslim.index <- support_general[treatment==3]
treatment_Cost.index <- support_general[treatment==4]
treatment_Sunni_Cost.index <- support_general[treatment==5]
treatment_Muslim_Cost.index <- support_general[treatment==6]


groups6 <- seq(1,6,1)

groups.index <- cbind(control.index, treatment_Sunni.index, treatment_Muslim.index, treatment_Cost.index, treatment_Sunni_Cost.index, treatment_Muslim_Cost.index)

means6 <- colMeans(groups.index, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses6 <- apply(groups.index, 2, se.mean)

uis6 <- means6 + 1.96*ses6
lis6 <- means6 - 1.96*ses6

data6 <- as.data.frame(cbind(condition, means6, ses6, uis6, lis6))


P6 <- ggplot(data6, aes(x=condition, y=means6)) +      
  geom_errorbar(aes(ymin=lis6, ymax=uis6),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(shape=21, size=3, fill = "black")+
  ylab("Mean index") +
  ggtitle("Figure 1.6. Composite Index") +
  coord_cartesian(ylim=c(-1, 2))+  
  #scale_y_continuous(breaks=seq(0,1, 0.1))+
  scale_x_discrete(limits = c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))+
  theme(axis.title.x = element_blank()) 



### CREATE MULTIPLE-PLOT FIGURE

multiplot(P1, P4, P2, P5, P3, P6, cols =3)

















################################################################
#### Hetereogeneous Effects ####################################
################################################################


######## BY CITY

control.stay.I <- shouldstay[treatment==1 & Gaziantep==0]
treatment_Sunni.stay.I <- shouldstay[treatment==2 & Gaziantep==0]
treatment_Muslim.stay.I <- shouldstay[treatment==3 & Gaziantep==0]
treatment_Cost.stay.I <- shouldstay[treatment==4 & Gaziantep==0]
treatment_Sunni_Cost.stay.I <- shouldstay[treatment==5 & Gaziantep==0]
treatment_Muslim_Cost.stay.I <- shouldstay[treatment==6 & Gaziantep==0]

groups <- seq(1,6,1)

groups.stayI <- cbind(control.stay.I, treatment_Sunni.stay.I, treatment_Muslim.stay.I, treatment_Cost.stay.I, treatment_Sunni_Cost.stay.I, treatment_Muslim_Cost.stay.I)

means1 <- colMeans(groups.stayI, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses1 <- apply(groups.stayI, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1


control.stay.G <- shouldstay[treatment==1 & Gaziantep==1]
treatment_Sunni.stay.G <- shouldstay[treatment==2 & Gaziantep==1]
treatment_Muslim.stay.G <- shouldstay[treatment==3 & Gaziantep==1]
treatment_Cost.stay.G <- shouldstay[treatment==4 & Gaziantep==1]
treatment_Sunni_Cost.stay.G <- shouldstay[treatment==5 & Gaziantep==1]
treatment_Muslim_Cost.stay.G <- shouldstay[treatment==6 & Gaziantep==1]

groups <- seq(1,6,1)

groups.stayG <- cbind(control.stay.G, treatment_Sunni.stay.G, treatment_Muslim.stay.G, treatment_Cost.stay.G, treatment_Sunni_Cost.stay.G, treatment_Muslim_Cost.stay.G)

means2 <- colMeans(groups.stayG, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.stayG, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


condition <- as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6))
levels(condition) <- c("Control","Sunni", "Muslim","Cost", 
                       "Sunni/Cost","Muslim/Cost")

group<-factor(c(1,2, 1,2, 1,2, 1,2, 1,2, 1,2))
levels(group) <- c("Istanbul","Gaziantep")

x1<-as.factor(1:12)

data1.C <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)

P1.C <- ggplot(data1.C, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Istanbul", "Gaziantep"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Istanbul", "Gaziantep"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 2.1. Refugees Can Stay\n") +
  coord_cartesian(ylim=c(0, 1))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))


### Figure 2.2

control.spend.I <- spend_ref[treatment==1 & Gaziantep==0]
treatment_Sunni.spend.I <- spend_ref[treatment==2 & Gaziantep==0]
treatment_Muslim.spend.I <- spend_ref[treatment==3 & Gaziantep==0]
treatment_Cost.spend.I <- spend_ref[treatment==4 & Gaziantep==0]
treatment_Sunni_Cost.spend.I <- spend_ref[treatment==5 & Gaziantep==0]
treatment_Muslim_Cost.spend.I <- spend_ref[treatment==6 & Gaziantep==0]

groups <- seq(1,6,1)

groups.spendI <- cbind(control.spend.I, treatment_Sunni.spend.I, treatment_Muslim.spend.I, treatment_Cost.spend.I, treatment_Sunni_Cost.spend.I, treatment_Muslim_Cost.spend.I)

means1 <- colMeans(groups.spendI, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses1 <- apply(groups.spendI, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1

control.spend.G <- spend_ref[treatment==1 & Gaziantep==1]
treatment_Sunni.spend.G <- spend_ref[treatment==2 & Gaziantep==1]
treatment_Muslim.spend.G <- spend_ref[treatment==3 & Gaziantep==1]
treatment_Cost.spend.G <- spend_ref[treatment==4 & Gaziantep==1]
treatment_Sunni_Cost.spend.G <- spend_ref[treatment==5 & Gaziantep==1]
treatment_Muslim_Cost.spend.G <- spend_ref[treatment==6 & Gaziantep==1]

groups <- seq(1,6,1)

groups.spendG <- cbind(control.spend.G, treatment_Sunni.spend.G, treatment_Muslim.spend.G, treatment_Cost.spend.G, treatment_Sunni_Cost.spend.G, treatment_Muslim_Cost.spend.G)

means2 <- colMeans(groups.spendG, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.spendG, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


condition <- as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6))
levels(condition) <- c("Control","Sunni", "Muslim","Cost", 
                       "Sunni/Cost","Muslim/Cost")

group<-factor(c(1,2, 1,2, 1,2, 1,2, 1,2, 1,2))
levels(group) <- c("Istanbul","Gaziantep")

x1<-as.factor(1:12)

data2.C <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)

P2.C <- ggplot(data2.C, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Istanbul", "Gaziantep"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Istanbul", "Gaziantep"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 2.2. Spend\n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))





### Figure 2.3

control.trust.I <- trust_ref[treatment==1 & Gaziantep==0]
treatment_Sunni.trust.I <- trust_ref[treatment==2 & Gaziantep==0]
treatment_Muslim.trust.I <- trust_ref[treatment==3 & Gaziantep==0]
treatment_Cost.trust.I <- trust_ref[treatment==4 & Gaziantep==0]
treatment_Sunni_Cost.trust.I<- trust_ref[treatment==5 & Gaziantep==0]
treatment_Muslim_Cost.trust.I <- trust_ref[treatment==6 & Gaziantep==0]

groups <- seq(1,6,1)

groups.trustI <- cbind(control.trust.I, treatment_Sunni.trust.I, treatment_Muslim.trust.I, treatment_Cost.trust.I, treatment_Sunni_Cost.trust.I, treatment_Muslim_Cost.trust.I)

means1 <- colMeans(groups.trustI, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses1 <- apply(groups.trustI, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1




control.trust.G <- trust_ref[treatment==1 & Gaziantep==1]
treatment_Sunni.trust.G <- trust_ref[treatment==2 & Gaziantep==1]
treatment_Muslim.trust.G <- trust_ref[treatment==3 & Gaziantep==1]
treatment_Cost.trust.G <- trust_ref[treatment==4 & Gaziantep==1]
treatment_Sunni_Cost.trust.G <- trust_ref[treatment==5 & Gaziantep==1]
treatment_Muslim_Cost.trust.G <- trust_ref[treatment==6 & Gaziantep==1]

groups <- seq(1,6,1)

groups.trustG <- cbind(control.trust.G, treatment_Sunni.trust.G, treatment_Muslim.trust.G, treatment_Cost.trust.G, treatment_Sunni_Cost.trust.G, treatment_Muslim_Cost.trust.G)

means2 <- colMeans(groups.trustG, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.trustG, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


condition <- as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6))
levels(condition) <- c("Control","Sunni", "Muslim","Cost", 
                       "Sunni/Cost","Muslim/Cost")

group<-factor(c(1,2, 1,2, 1,2, 1,2, 1,2, 1,2))
levels(group) <- c("Istanbul","Gaziantep")

x1<-as.factor(1:12)

data3.C <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)

P3.C <- ggplot(data3.C, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Istanbul", "Gaziantep"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Istanbul", "Gaziantep"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 2.3. Trust\n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



### Figure 2.4

control.n.I <- neighborhood_ref[treatment==1 & Gaziantep==0]
treatment_Sunni.n.I <- neighborhood_ref[treatment==2 & Gaziantep==0]
treatment_Muslim.n.I <- neighborhood_ref[treatment==3 & Gaziantep==0]
treatment_Cost.n.I <- neighborhood_ref[treatment==4 & Gaziantep==0]
treatment_Sunni_Cost.n.I <- neighborhood_ref[treatment==5 & Gaziantep==0]
treatment_Muslim_Cost.n.I <- neighborhood_ref[treatment==6 & Gaziantep==0]

groups <- seq(1,6,1)

groups.nI <- cbind(control.n.I, treatment_Sunni.n.I, treatment_Muslim.n.I, treatment_Cost.n.I, treatment_Sunni_Cost.n.I, treatment_Muslim_Cost.n.I)

means1 <- colMeans(groups.nI, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses1 <- apply(groups.nI, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1



control.n.G <- neighborhood_ref[treatment==1 & Gaziantep==1]
treatment_Sunni.n.G <- neighborhood_ref[treatment==2 & Gaziantep==1]
treatment_Muslim.n.G <- neighborhood_ref[treatment==3 & Gaziantep==1]
treatment_Cost.n.G <- neighborhood_ref[treatment==4 & Gaziantep==1]
treatment_Sunni_Cost.n.G <- neighborhood_ref[treatment==5 & Gaziantep==1]
treatment_Muslim_Cost.n.G <- neighborhood_ref[treatment==6 & Gaziantep==1]

groups <- seq(1,6,1)

groups.nG <- cbind(control.n.G, treatment_Sunni.n.G, treatment_Muslim.n.G, treatment_Cost.n.G, treatment_Sunni_Cost.n.G, treatment_Muslim_Cost.n.G)

means2 <- colMeans(groups.nG, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.nG, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


condition <- as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6))
levels(condition) <- c("Control","Sunni", "Muslim","Cost", 
                       "Sunni/Cost","Muslim/Cost")

group<-factor(c(1,2, 1,2, 1,2, 1,2, 1,2, 1,2))
levels(group) <- c("Istanbul","Gaziantep")

x1<-as.factor(1:12)

data4.C <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)

P4.C <- ggplot(data4.C, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Istanbul", "Gaziantep"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Istanbul", "Gaziantep"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 2.4. Acceptance of the refugees \n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



### Figure 2.5

control.d.I <- donation[treatment==1 & Gaziantep==0]
treatment_Sunni.d.I <- donation[treatment==2 & Gaziantep==0]
treatment_Muslim.d.I <- donation[treatment==3 & Gaziantep==0]
treatment_Cost.d.I <- donation[treatment==4 & Gaziantep==0]
treatment_Sunni_Cost.d.I <- donation[treatment==5 & Gaziantep==0]
treatment_Muslim_Cost.d.I <- donation[treatment==6 & Gaziantep==0]

groups <- seq(1,6,1)

groups.dI <- cbind(control.d.I, treatment_Sunni.d.I, treatment_Muslim.d.I, treatment_Cost.d.I, treatment_Sunni_Cost.d.I, treatment_Muslim_Cost.d.I)

means1 <- colMeans(groups.dI, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses1 <- apply(groups.dI, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1


control.d.G <- donation[treatment==1 & Gaziantep==1]
treatment_Sunni.d.G <- donation[treatment==2 & Gaziantep==1]
treatment_Muslim.d.G <- donation[treatment==3 & Gaziantep==1]
treatment_Cost.d.G <- donation[treatment==4 & Gaziantep==1]
treatment_Sunni_Cost.d.G <- donation[treatment==5 & Gaziantep==1]
treatment_Muslim_Cost.d.G <- donation[treatment==6 & Gaziantep==1]

groups <- seq(1,6,1)

groups.dG <- cbind(control.d.G, treatment_Sunni.d.G, treatment_Muslim.d.G, treatment_Cost.d.G, treatment_Sunni_Cost.d.G, treatment_Muslim_Cost.d.G)

means2 <- colMeans(groups.dG, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.dG, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


condition <- as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6))
levels(condition) <- c("Control","Sunni", "Muslim","Cost", 
                       "Sunni/Cost","Muslim/Cost")

group<-factor(c(1,2, 1,2, 1,2, 1,2, 1,2, 1,2))
levels(group) <- c("Istanbul","Gaziantep")

x1<-as.factor(1:12)

data5.C <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)

P5.C <- ggplot(data5.C, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Istanbul", "Gaziantep"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Istanbul", "Gaziantep"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 2.5. Donation\n") +
  coord_cartesian(ylim=c(0, 4))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



## Figure 2.6

control.index.I <- support_general[treatment==1 & Gaziantep==0]
treatment_Sunni.index.I <- support_general[treatment==2 & Gaziantep==0]
treatment_Muslim.index.I <-support_general[treatment==3 & Gaziantep==0]
treatment_Cost.index.I <- support_general[treatment==4 & Gaziantep==0]
treatment_Sunni_Cost.index.I <- support_general[treatment==5 & Gaziantep==0]
treatment_Muslim_Cost.index.I <- support_general[treatment==6 & Gaziantep==0]

groups <- seq(1,6,1)

groups.indexI <- cbind(control.index.I, treatment_Sunni.index.I, treatment_Muslim.index.I, treatment_Cost.index.I, treatment_Sunni_Cost.index.I, treatment_Muslim_Cost.index.I)

means1 <- colMeans(groups.indexI, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses1 <- apply(groups.indexI, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1


control.index.G <- support_general[treatment==1 & Gaziantep==1]
treatment_Sunni.index.G <- support_general[treatment==2 & Gaziantep==1]
treatment_Muslim.index.G <- support_general[treatment==3 & Gaziantep==1]
treatment_Cost.index.G <- support_general[treatment==4 & Gaziantep==1]
treatment_Sunni_Cost.index.G <- support_general[treatment==5 & Gaziantep==1]
treatment_Muslim_Cost.index.G <- support_general[treatment==6 & Gaziantep==1]

groups <- seq(1,6,1)

groups.indexG <- cbind(control.index.G, treatment_Sunni.index.G, treatment_Muslim.index.G, treatment_Cost.index.G, treatment_Sunni_Cost.index.G, treatment_Muslim_Cost.index.G)

means2 <- colMeans(groups.indexG, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.indexG, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


condition <- as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6))
levels(condition) <- c("Control","Sunni", "Muslim","Cost", 
                       "Sunni/Cost","Muslim/Cost")

group<-factor(c(1,2, 1,2, 1,2, 1,2, 1,2, 1,2))
levels(group) <- c("Istanbul","Gaziantep")

x1<-as.factor(1:12)

data6.C <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)

P6.C <- ggplot(data6.C, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Istanbul", "Gaziantep"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Istanbul", "Gaziantep"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 2.6. Composite Index\n") +
  coord_cartesian(ylim=c(-1, 1))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))

multiplot(P1.C, P4.C, P2.C, P5.C, P3.C, P6.C, cols =3) 





######## BY CONTACT


### Figure 3.1

control.stay.c <- shouldstay[treatment==1 & contact_dummy==1]
treatment_Sunni.stay.c <- shouldstay[treatment==2 & contact_dummy==1]
treatment_Muslim.stay.c <- shouldstay[treatment==3 & contact_dummy==1]
treatment_Cost.stay.c <- shouldstay[treatment==4 & contact_dummy==1]
treatment_Sunni_Cost.stay.c <- shouldstay[treatment==5 & contact_dummy==1]
treatment_Muslim_Cost.stay.c <- shouldstay[treatment==6 & contact_dummy==1]

groups <- seq(1,6,1)

groups.stay.c <- cbind(control.stay.c, treatment_Sunni.stay.c, treatment_Muslim.stay.c, treatment_Cost.stay.c, treatment_Sunni_Cost.stay.c, treatment_Muslim_Cost.stay.c)

means3.1_1 <- colMeans(groups.stay.c, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses3.1_1 <- apply(groups.stay.c, 2, se.mean)

uis3.1_1 <- means3.1_1 + 1.96*ses3.1_1
lis3.1_1 <- means3.1_1 - 1.96*ses3.1_1


control.stay.nc <- shouldstay[treatment==1 & contact_dummy==0]
treatment_Sunni.stay.nc <- shouldstay[treatment==2 & contact_dummy==0]
treatment_Muslim.stay.nc <- shouldstay[treatment==3 & contact_dummy==0]
treatment_Cost.stay.nc <- shouldstay[treatment==4 & contact_dummy==0]
treatment_Sunni_Cost.stay.nc <- shouldstay[treatment==5 & contact_dummy==0]
treatment_Muslim_Cost.stay.nc <- shouldstay[treatment==6 & contact_dummy==0]


groups <- seq(1,6,1)

groups.stay.nc <- cbind(control.stay.nc, treatment_Sunni.stay.nc, treatment_Muslim.stay.nc, treatment_Cost.stay.nc, treatment_Sunni_Cost.stay.nc, treatment_Muslim_Cost.stay.nc)

means3.1_2 <- colMeans(groups.stay.nc, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses3.1_2 <- apply(groups.stay.nc, 2, se.mean)

uis3.1_2 <- means3.1_2  + 1.96*ses3.1_2 
lis3.1_2 <- means3.1_2  - 1.96*ses3.1_2 


means3.1 <- c(means3.1_1, means3.1_2)
means3.1 <- means3.1[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses3.1 <- c(ses3.1_1, ses3.1_2)
ses3.1 <- ses3.1[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis3.1 <- c(uis3.1_1, uis3.1_2)
uis3.1 <- uis3.1[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis3.1 <- c(lis3.1_1, lis3.1_2)
lis3.1 <- lis3.1[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


condition <- as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6))
levels(condition) <- c("Control","Sunni", "Muslim","Cost", 
                       "Sunni/Cost","Muslim/Cost")


group<-factor(c(1,2, 1,2, 1,2, 1,2, 1,2, 1,2))
levels(group) <- c("Contact","No Contact")

x1<-as.factor(1:12)

data1.c <- as.data.frame(cbind(x1,condition, means3.1, ses3.1, uis3.1, lis3.1, group),
                       use.value.labels= T)

P1.c <- ggplot(data1.c, aes(x=x1, y=means3.1,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis3.1, ymax=uis3.1, width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Contact", "No Contact"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Contact", "No Contact"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 3.1. Refugees Can Stay\n") +
  coord_cartesian(ylim=c(0, 1))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))


### Figure 3.2

control.spend.c <- spend_ref[treatment==1 & contact_dummy==1]
treatment_Sunni.spend.c <- spend_ref[treatment==2 & contact_dummy==1]
treatment_Muslim.spend.c <- spend_ref[treatment==3 & contact_dummy==1]
treatment_Cost.spend.c <- spend_ref[treatment==4 & contact_dummy==1]
treatment_Sunni_Cost.spend.c <- spend_ref[treatment==5 & contact_dummy==1]
treatment_Muslim_Cost.spend.c <- spend_ref[treatment==6 & contact_dummy==1]


groups <- seq(1,6,1)

groups.spend.c <- cbind(control.spend.c, treatment_Sunni.spend.c, treatment_Muslim.spend.c, treatment_Cost.spend.c, treatment_Sunni_Cost.spend.c, treatment_Muslim_Cost.spend.c)

means3.2_1 <- colMeans(groups.spend.c, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3.2_1 <- apply(groups.spend.c, 2, se.mean)

uis3.2_1 <- means3.2_1 + 1.96*ses3.2_1
lis3.2_1 <- means3.2_1 - 1.96*ses3.2_1


control.spend.nc <- spend_ref[treatment==1 & contact_dummy==0]
treatment_Sunni.spend.nc <- spend_ref[treatment==2 & contact_dummy==0]
treatment_Muslim.spend.nc <- spend_ref[treatment==3 & contact_dummy==0]
treatment_Cost.spend.nc <- spend_ref[treatment==4 & contact_dummy==0]
treatment_Sunni_Cost.spend.nc <- spend_ref[treatment==5 & contact_dummy==0]
treatment_Muslim_Cost.spend.nc <- spend_ref[treatment==6 & contact_dummy==0]


groups <- seq(1,6,1)

groups.spend.nc <- cbind(control.spend.nc, treatment_Sunni.spend.nc, treatment_Muslim.spend.nc, treatment_Cost.spend.nc, treatment_Sunni_Cost.spend.nc, treatment_Muslim_Cost.spend.nc)

means3.2_2 <- colMeans(groups.spend.nc, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses3.2_2 <- apply(groups.spend.nc, 2, se.mean)

uis3.2_2 <- means3.2_2 + 1.96*ses3.2_2
lis3.2_2 <- means3.2_2 - 1.96*ses3.2_2



means3.2 <- c(means3.2_1, means3.2_2)
means3.2 <- means3.2[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses3.2 <- c(ses3.2_1, ses3.2_2)
ses3.2 <- ses3.2[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis3.2 <- c(uis3.2_1, uis3.2_2)
uis3.2 <- uis3.2[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis3.2 <- c(lis3.2_1, lis3.2_2)
lis3.2 <- lis3.2[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]




data2.c <- as.data.frame(cbind(x1,condition, means3.2, ses3.2, uis3.2, lis3.2, group),
                         use.value.labels= T)



P2.c <- ggplot(data2.c, aes(x=x1, y=means3.2,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis3.2, ymax=uis3.2,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Contact", "No Contact"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Contact", "No Contact"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 3.2. Support for Government Spending\n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))





### Figure 3.3

control.trust.c <- trust_ref[treatment==1 & contact_dummy==1]
treatment_Sunni.trust.c <- trust_ref[treatment==2 & contact_dummy==1]
treatment_Muslim.trust.c <- trust_ref[treatment==3 & contact_dummy==1]
treatment_Cost.trust.c <- trust_ref[treatment==4 & contact_dummy==1]
treatment_Sunni_Cost.trust.c <- trust_ref[treatment==5 & contact_dummy==1]
treatment_Muslim_Cost.trust.c <- trust_ref[treatment==6 & contact_dummy==1]

groups <- seq(1,6,1)

groups.trust.c <- cbind(control.trust.c, treatment_Sunni.trust.c, treatment_Muslim.trust.c, treatment_Cost.trust.c, treatment_Sunni_Cost.trust.c, treatment_Muslim_Cost.trust.c)

means3.3_1 <- colMeans(groups.trust.c, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3.3_1 <- apply(groups.trust.c, 2, se.mean)

uis3.3_1 <- means3.3_1 + 1.96*ses3.3_1
lis3.3_1 <- means3.3_1 - 1.96*ses3.3_1



control.trust.nc <- trust_ref[treatment==1 & contact_dummy==0]
treatment_Sunni.trust.nc <- trust_ref[treatment==2 & contact_dummy==0]
treatment_Muslim.trust.nc <- trust_ref[treatment==3 & contact_dummy==0]
treatment_Cost.trust.nc <- trust_ref[treatment==4 & contact_dummy==0]
treatment_Sunni_Cost.trust.nc <- trust_ref[treatment==5 & contact_dummy==0]
treatment_Muslim_Cost.trust.nc <- trust_ref[treatment==6 & contact_dummy==0]

groups <- seq(1,6,1)

groups.trust.nc <- cbind(control.trust.nc, treatment_Sunni.trust.nc, treatment_Muslim.trust.nc, treatment_Cost.trust.nc, treatment_Sunni_Cost.trust.nc, treatment_Muslim_Cost.trust.nc)

means3.3_2 <- colMeans(groups.trust.nc, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses3.3_2 <- apply(groups.trust.nc, 2, se.mean)

uis3.3_2 <- means3.3_2 + 1.96*ses3.3_2
lis3.3_2 <- means3.3_2 - 1.96*ses3.3_2


means3.3 <- c(means3.3_1, means3.3_2)
means3.3 <- means3.3[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses3.3 <- c(ses3.3_1, ses3.3_2)
ses3.3 <- ses3.3[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis3.3 <- c(uis3.3_1, uis3.3_2)
uis3.3 <- uis3.3[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis3.3 <- c(lis3.3_1, lis3.3_2)
lis3.3 <- lis3.3[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


data3.c <- as.data.frame(cbind(x1,condition, means3.3, ses3.3, uis3.3, lis3.3, group),use.value.labels= T)


P3.c <- ggplot(data3.c, aes(x=x1, y=means3.3,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis3.3, ymax=uis3.3,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Contact", "No Contact"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Contact", "No Contact"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 3.3. Trust\n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



### Figure 3.4

control.n.c <- neighborhood_ref[treatment==1 & contact_dummy==1]
treatment_Sunni.n.c <- neighborhood_ref[treatment==2 & contact_dummy==1]
treatment_Muslim.n.c <- neighborhood_ref[treatment==3 & contact_dummy==1]
treatment_Cost.n.c <- neighborhood_ref[treatment==4 & contact_dummy==1]
treatment_Sunni_Cost.n.c <- neighborhood_ref[treatment==5 & contact_dummy==1]
treatment_Muslim_Cost.n.c <- neighborhood_ref[treatment==6 & contact_dummy==1]

groups <- seq(1,6,1)

groups.n.c <- cbind(control.n.c, treatment_Sunni.n.c, treatment_Muslim.n.c, treatment_Cost.n.c, treatment_Sunni_Cost.n.c, treatment_Muslim_Cost.n.c)

means3.4_1 <- colMeans(groups.n.c, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3.4_1 <- apply(groups.n.c, 2, se.mean)

uis3.4_1 <- means3.4_1 + 1.96*ses3.4_1
lis3.4_1 <- means3.4_1 - 1.96*ses3.4_1




control.n.nc <- neighborhood_ref[treatment==1 & contact_dummy==0]
treatment_Sunni.n.nc <- neighborhood_ref[treatment==2 & contact_dummy==0]
treatment_Muslim.n.nc <- neighborhood_ref[treatment==3 & contact_dummy==0]
treatment_Cost.n.nc <- neighborhood_ref[treatment==4 & contact_dummy==0]
treatment_Sunni_Cost.n.nc <- neighborhood_ref[treatment==5 & contact_dummy==0]
treatment_Muslim_Cost.n.nc <- neighborhood_ref[treatment==6 & contact_dummy==0]


groups <- seq(1,6,1)

groups.n.nc <- cbind(control.n.nc, treatment_Sunni.n.nc, treatment_Muslim.n.nc, treatment_Cost.n.nc, treatment_Sunni_Cost.n.nc, treatment_Muslim_Cost.n.nc)

means3.4_2 <- colMeans(groups.n.nc, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3.4_2 <- apply(groups.n.nc, 2, se.mean)

uis3.4_2 <- means3.4_2 + 1.96*ses3.4_2
lis3.4_2 <- means3.4_2 - 1.96*ses3.4_2



means3.4 <- c(means3.4_1,means3.4_2)
means3.4 <- means3.4[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses3.4 <- c(ses3.4_1, ses3.4_2)
ses3.4 <- ses3.4[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis3.4 <- c(uis3.4_1, uis3.4_2)
uis3.4 <- uis3.4[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis3.4 <- c(lis3.4_1, lis3.4_2)
lis3.4 <- lis3.4[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]



data4.c <- as.data.frame(cbind(x1,condition, means3.4, ses3.4, uis3.4, lis3.4, group),
                         use.value.labels= T)



P4.c <- ggplot(data4.c, aes(x=x1, y=means3.4,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis3.4, ymax=uis3.4,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Contact", "No Contact"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Contact", "No Contact"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 3.4. Acceptance of the refugees \n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



### Figure 3.5

control.don.c <- donation[treatment==1 & contact_dummy==1]
treatment_Sunni.don.c <- donation[treatment==2 & contact_dummy==1]
treatment_Muslim.don.c <- donation[treatment==3 & contact_dummy==1]
treatment_Cost.don.c <- donation[treatment==4 & contact_dummy==1]
treatment_Sunni_Cost.don.c <- donation[treatment==5 & contact_dummy==1]
treatment_Muslim_Cost.don.c <- donation[treatment==6 & contact_dummy==1]

groups <- seq(1,6,1)

groups.don.c <- cbind(control.don.c, treatment_Sunni.don.c, treatment_Muslim.don.c, treatment_Cost.don.c, treatment_Sunni_Cost.don.c, treatment_Muslim_Cost.don.c)

means3.5_1 <- colMeans(groups.don.c, na.rm=TRUE)


se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3.5_1 <- apply(groups.don.c, 2, se.mean)

uis3.5_1 <- means3.5_1 + 1.96*ses3.5_1
lis3.5_1 <- means3.5_1 - 1.96*ses3.5_1



control.don.nc <- donation[treatment==1 & contact_dummy==0]
treatment_Sunni.don.nc <- donation[treatment==2 & contact_dummy==0]
treatment_Muslim.don.nc <- donation[treatment==3 & contact_dummy==0]
treatment_Cost.don.nc <- donation[treatment==4 & contact_dummy==0]
treatment_Sunni_Cost.don.nc <- donation[treatment==5 & contact_dummy==0]
treatment_Muslim_Cost.don.nc <- donation[treatment==6 & contact_dummy==0]

groups <- seq(1,6,1)

groups.don.nc <- cbind(control.don.nc, treatment_Sunni.don.nc, treatment_Muslim.don.nc, treatment_Cost.don.nc, treatment_Sunni_Cost.don.nc, treatment_Muslim_Cost.don.nc)

means3.5_2 <- colMeans(groups.don.nc, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3.5_2 <- apply(groups.don.nc, 2, se.mean)

uis3.5_2 <- means3.5_2 + 1.96*ses3.5_2
lis3.5_2 <- means3.5_2 - 1.96*ses3.5_2



means3.5 <- c(means3.5_1, means3.5_2)
means3.5 <- means3.5[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses3.5 <- c(ses3.5_1, ses3.5_2)
ses3.5 <- ses3.5[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis3.5 <- c(uis3.5_1, uis3.5_2)
uis3.5 <- uis3.5[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis3.5 <- c(lis3.5_1, lis3.5_2)
lis3.5 <- lis3.5[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]




data5.c <- as.data.frame(cbind(x1,condition, means3.5, ses3.5, uis3.5, lis3.5, group),
                         use.value.labels= T)


P5.c <- ggplot(data5.c, aes(x=x1, y=means3.5,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis3.5, ymax=uis3.5, width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Contact", "No contact"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Contact", "No contact"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 3.5. Donation\n") +
  coord_cartesian(ylim=c(0, 4))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



## Figure 3.6

control.index.c <- support_general[treatment==1 & contact_dummy==1]
treatment_Sunni.index.c <- support_general[treatment==2 & contact_dummy==1]
treatment_Muslim.index.c <- support_general[treatment==3 & contact_dummy==1]
treatment_Cost.index.c <- support_general[treatment==4 & contact_dummy==1]
treatment_Sunni_Cost.index.c <- support_general[treatment==5 & contact_dummy==1]
treatment_Muslim_Cost.index.c <- support_general[treatment==6 & contact_dummy==1]

groups <- seq(1,6,1)

groups.index.c <- cbind(control.index.c, treatment_Sunni.index.c, treatment_Muslim.index.c, treatment_Cost.index.c, treatment_Sunni_Cost.index.c, treatment_Muslim_Cost.index.c)

means3.6_1 <- colMeans(groups.index.c, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3.6_1 <- apply(groups.index.c, 2, se.mean)

uis3.6_1 <- means3.6_1 + 1.96*ses3.6_1
lis3.6_1 <- means3.6_1 - 1.96*ses3.6_1



control.index.nc <- support_general[treatment==1 & contact_dummy==0]
treatment_Sunni.index.nc <- support_general[treatment==2 & contact_dummy==0]
treatment_Muslim.index.nc <- support_general[treatment==3 & contact_dummy==0]
treatment_Cost.index.nc <- support_general[treatment==4 & contact_dummy==0]
treatment_Sunni_Cost.index.nc <- support_general[treatment==5 & contact_dummy==0]
treatment_Muslim_Cost.index.nc <- support_general[treatment==6 & contact_dummy==0]

groups <- seq(1,6,1)

groups.index.nc <- cbind( control.index.nc, treatment_Sunni.index.nc, treatment_Muslim.index.nc, treatment_Cost.index.nc, treatment_Sunni_Cost.index.nc, treatment_Muslim_Cost.index.nc)

means3.6_2 <- colMeans(groups.index.nc, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses3.6_2  <- apply(groups.index.nc, 2, se.mean)

uis3.6_2  <- means3.6_2  + 1.96*ses3.6_2 
lis3.6_2  <- means3.6_2  - 1.96*ses3.6_2 


means3.6 <- c(means3.6_1,means3.6_2 )
means3.6 <- means3.6[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses3.6 <- c(ses3.6_1 , ses3.6_2 )
ses3.6 <- ses3.6[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis3.6 <- c(uis3.6_1 , uis3.6_2 )
uis3.6 <- uis3.6[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis3.6 <- c(lis3.6_1, lis3.6_2 )
lis3.6 <- lis3.6[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


data6.c <- as.data.frame(cbind(x1,condition, means3.6, ses3.6, uis3.6, lis3.6, group),
                         use.value.labels= T)


P6.c <- ggplot(data6.c, aes(x=x1, y=means3.6,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis3.6, ymax=uis3.6,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Contact", "No Contact"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Contact", "No Contact"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 3.6. Composite Index\n") +
  coord_cartesian(ylim=c(-1, 1))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))

multiplot(P1.c, P4.c, P2.c, P5.c, P3.c, P6.c, cols =3) 

















######## BY RELIGIOSITY


### Figure 4.1

control.stay.r <- shouldstay[treatment==1 & religiosity_dummy==1]
treatment_Sunni.stay.r <- shouldstay[treatment==2 & religiosity_dummy==1]
treatment_Muslim.stay.r <- shouldstay[treatment==3 & religiosity_dummy==1]
treatment_Cost.stay.r <- shouldstay[treatment==4 & religiosity_dummy==1]
treatment_Sunni_Cost.stay.r <- shouldstay[treatment==5 & religiosity_dummy==1]
treatment_Muslim_Cost.stay.r <- shouldstay[treatment==6 & religiosity_dummy==1]

groups <- seq(1,6,1)

groups.stay.r <- cbind(control.stay.r, treatment_Sunni.stay.r, treatment_Muslim.stay.r, treatment_Cost.stay.r, treatment_Sunni_Cost.stay.r, treatment_Muslim_Cost.stay.r)

means1 <- colMeans(groups.stay.r, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses1 <- apply(groups.stay.r, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1





control.stay.nr <- shouldstay[treatment==1 & religiosity_dummy==0]
treatment_Sunni.stay.nr <- shouldstay[treatment==2 & religiosity_dummy==0]
treatment_Muslim.stay.nr <- shouldstay[treatment==3 & religiosity_dummy==0]
treatment_Cost.stay.nr <- shouldstay[treatment==4 & religiosity_dummy==0]
treatment_Sunni_Cost.stay.nr <- shouldstay[treatment==5 & religiosity_dummy==0]
treatment_Muslim_Cost.stay.nr <- shouldstay[treatment==6 & religiosity_dummy==0]


groups <- seq(1,6,1)

groups.stay.nr <- cbind(control.stay.nr, treatment_Sunni.stay.nr, treatment_Muslim.stay.nr, treatment_Cost.stay.nr, treatment_Sunni_Cost.stay.nr, treatment_Muslim_Cost.stay.nr)

means2 <- colMeans(groups.stay.nr, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.stay.nr, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


condition <- as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6))
levels(condition) <- c("Control","Sunni", "Muslim","Cost", 
                       "Sunni/Cost","Muslim/Cost")


group<-factor(c(1,2, 1,2, 1,2, 1,2, 1,2, 1,2))
levels(group) <- c("Highly Religious ","Less Religious")

x1<-as.factor(1:12)

data1.r <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)

P1.r <- ggplot(data1.r, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Highly Religious", "Less Religious"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Highly Religious", "Less Religious"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 4.1. Refugees Can Stay\n") +
  coord_cartesian(ylim=c(0, 1))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))


### Figure 4.2

control.spend.r <- spend_ref[treatment==1 & religiosity_dummy==1]
treatment_Sunni.spend.r <- spend_ref[treatment==2 & religiosity_dummy==1]
treatment_Muslim.spend.r <- spend_ref[treatment==3 & religiosity_dummy==1]
treatment_Cost.spend.r <- spend_ref[treatment==4 & religiosity_dummy==1]
treatment_Sunni_Cost.spend.r <- spend_ref[treatment==5 & religiosity_dummy==1]
treatment_Muslim_Cost.spend.r <- spend_ref[treatment==6 & religiosity_dummy==1]


groups <- seq(1,6,1)

groups.spend.r <- cbind(control.spend.r, treatment_Sunni.spend.r, treatment_Muslim.spend.r, treatment_Cost.spend.r, treatment_Sunni_Cost.spend.r, treatment_Muslim_Cost.spend.r)

means1 <- colMeans(groups.spend.r, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses1 <- apply(groups.spend.r, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1


control.spend.nr <- spend_ref[treatment==1 & religiosity_dummy==0]
treatment_Sunni.spend.nr <- spend_ref[treatment==2 & religiosity_dummy==0]
treatment_Muslim.spend.nr <- spend_ref[treatment==3 & religiosity_dummy==0]
treatment_Cost.spend.nr <- spend_ref[treatment==4 & religiosity_dummy==0]
treatment_Sunni_Cost.spend.nr <- spend_ref[treatment==5 & religiosity_dummy==0]
treatment_Muslim_Cost.spend.nr <- spend_ref[treatment==6 & religiosity_dummy==0]


groups <- seq(1,6,1)

groups.spend.nr <- cbind(control.spend.nr, treatment_Sunni.spend.nr, treatment_Muslim.spend.nr, treatment_Cost.spend.nr, treatment_Sunni_Cost.spend.nr, treatment_Muslim_Cost.spend.nr)

means2 <- colMeans(groups.spend.nr, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.spend.nr, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2



means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]




data2.r <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)


P2.r <- ggplot(data2.r, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Highly Religious", "Less Religious"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Highly Religious", "Less Religious"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 4.2. Support for Government Spending\n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



### Figure 4.3

control.trust.r <- trust_ref[treatment==1 & religiosity_dummy==0]
treatment_Sunni.trust.r <- trust_ref[treatment==2 & religiosity_dummy==0]
treatment_Muslim.trust.r <- trust_ref[treatment==3 & religiosity_dummy==0]
treatment_Cost.trust.r <- trust_ref[treatment==4 & religiosity_dummy==0]
treatment_Sunni_Cost.trust.r <- trust_ref[treatment==5 & religiosity_dummy==0]
treatment_Muslim_Cost.trust.r <- trust_ref[treatment==6 & religiosity_dummy==0]

groups <- seq(1,6,1)

groups.trust.r <- cbind(control.trust.r, treatment_Sunni.trust.r, treatment_Muslim.trust.r, treatment_Cost.trust.r, treatment_Sunni_Cost.trust.r, treatment_Muslim_Cost.trust.r)

means1 <- colMeans(groups.trust.r, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses1 <- apply(groups.trust.r, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1



control.trust.nr <- trust_ref[treatment==1 & religiosity_dummy==0]
treatment_Sunni.trust.nr <- trust_ref[treatment==2 & religiosity_dummy==0]
treatment_Muslim.trust.nr <- trust_ref[treatment==3 & religiosity_dummy==0]
treatment_Cost.trust.nr <- trust_ref[treatment==4 & religiosity_dummy==0]
treatment_Sunni_Cost.trust.nr <- trust_ref[treatment==5 & religiosity_dummy==0]
treatment_Muslim_Cost.trust.nr <- trust_ref[treatment==6 & religiosity_dummy==0]

groups <- seq(1,6,1)

groups.trust.nr <- cbind(control.trust.nr, treatment_Sunni.trust.nr, treatment_Muslim.trust.nr, treatment_Cost.trust.nr, treatment_Sunni_Cost.trust.nr, treatment_Muslim_Cost.trust.nr)

means2 <- colMeans(groups.trust.nr, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}

ses2 <- apply(groups.trust.nr, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


data3.r <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),use.value.labels= T)


P3.r <- ggplot(data3.r, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Highly Religious", "Less Religious"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Highly Religious", "Less Religious"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 4.3. Trust\n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



### Figure 4.4

control.n.r <- neighborhood_ref[treatment==1 & religiosity_dummy==1]
treatment_Sunni.n.r <- neighborhood_ref[treatment==2 & religiosity_dummy==1]
treatment_Muslim.n.r <- neighborhood_ref[treatment==3 & religiosity_dummy==1]
treatment_Cost.n.r <- neighborhood_ref[treatment==4 & religiosity_dummy==1]
treatment_Sunni_Cost.n.r <- neighborhood_ref[treatment==5 & religiosity_dummy==1]
treatment_Muslim_Cost.n.r <- neighborhood_ref[treatment==6 & religiosity_dummy==1]

groups <- seq(1,6,1)

groups.n.r <- cbind(control.n.r, treatment_Sunni.n.r, treatment_Muslim.n.r, treatment_Cost.n.r, treatment_Sunni_Cost.n.r, treatment_Muslim_Cost.n.r)

means1 <- colMeans(groups.n.r, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}



ses1 <- apply(groups.n.r, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1




control.n.nr <- neighborhood_ref[treatment==1 & religiosity_dummy==0]
treatment_Sunni.n.nr <- neighborhood_ref[treatment==2 & religiosity_dummy==0]
treatment_Muslim.n.nr <- neighborhood_ref[treatment==3 & religiosity_dummy==0]
treatment_Cost.n.nr <- neighborhood_ref[treatment==4 & religiosity_dummy==0]
treatment_Sunni_Cost.n.nr <- neighborhood_ref[treatment==5 & religiosity_dummy==0]
treatment_Muslim_Cost.n.nr <- neighborhood_ref[treatment==6 & religiosity_dummy==0]


groups <- seq(1,6,1)

groups.n.nr <- cbind(control.n.nr, treatment_Sunni.n.nr, treatment_Muslim.n.nr, treatment_Cost.n.nr, treatment_Sunni_Cost.n.nr, treatment_Muslim_Cost.n.nr)

means2 <- colMeans(groups.n.nr, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses2 <- apply(groups.n.nr, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2



means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]



data4.r <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)



P4.r <- ggplot(data4.r, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Highly Religious", "Less Religious"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Highly Religious", "Less Religious"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 4.4. Acceptance of the refugees \n") +
  coord_cartesian(ylim=c(1, 3))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))





### Figure 4.5

control.don.r <- donation[treatment==1 & religiosity_dummy==1]
treatment_Sunni.don.r <- donation[treatment==2 & religiosity_dummy==1]
treatment_Muslim.don.r <- donation[treatment==3 & religiosity_dummy==1]
treatment_Cost.don.r <- donation[treatment==4 & religiosity_dummy==1]
treatment_Sunni_Cost.don.r <- donation[treatment==5 & religiosity_dummy==1]
treatment_Muslim_Cost.don.r <- donation[treatment==6 & religiosity_dummy==1]

groups <- seq(1,6,1)

groups.don.r <- cbind(control.don.r, treatment_Sunni.don.r, treatment_Muslim.don.r, treatment_Cost.don.r, treatment_Sunni_Cost.don.r, treatment_Muslim_Cost.don.r)

means1 <- colMeans(groups.don.r, na.rm=TRUE)


se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses1 <- apply(groups.don.r, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1


control.don.nr <- donation[treatment==1 & religiosity_dummy==0]
treatment_Sunni.don.nr <- donation[treatment==2 & religiosity_dummy==0]
treatment_Muslim.don.nr <- donation[treatment==3 & religiosity_dummy==0]
treatment_Cost.don.nr <- donation[treatment==4 & religiosity_dummy==0]
treatment_Sunni_Cost.don.nr <- donation[treatment==5 & religiosity_dummy==0]
treatment_Muslim_Cost.don.nr <- donation[treatment==6 & religiosity_dummy==0]

groups <- seq(1,6,1)

groups.don.nr <- cbind(control.don.nr, treatment_Sunni.don.nr, treatment_Muslim.don.nr, treatment_Cost.don.nr, treatment_Sunni_Cost.don.nr, treatment_Muslim_Cost.don.nr)

means2 <- colMeans(groups.don.nr, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses2 <- apply(groups.don.nr, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2

means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]




data5.r <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)


P5.r <- ggplot(data5.r, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Highly Religious", "Less Religious"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Highly Religious", "Less Religious"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 4.5. Donation\n") +
  coord_cartesian(ylim=c(0, 4))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))



## Figure 3.6

control.index.r <- support_general[treatment==1 & religiosity_dummy==1]
treatment_Sunni.index.r <- support_general[treatment==2 & religiosity_dummy==1]
treatment_Muslim.index.r <- support_general[treatment==3 & religiosity_dummy==1]
treatment_Cost.index.r <- support_general[treatment==4 & religiosity_dummy==1]
treatment_Sunni_Cost.index.r <- support_general[treatment==5 & religiosity_dummy==1]
treatment_Muslim_Cost.index.r <- support_general[treatment==6 & religiosity_dummy==1]

groups <- seq(1,6,1)

groups.index.r <- cbind(control.index.r, treatment_Sunni.index.r, treatment_Muslim.index.r, treatment_Cost.index.r, treatment_Sunni_Cost.index.r, treatment_Muslim_Cost.index.r)

means1 <- colMeans(groups.index.r, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses1 <- apply(groups.index.r, 2, se.mean)

uis1 <- means1 + 1.96*ses1
lis1 <- means1 - 1.96*ses1



control.index.nr <- support_general[treatment==1 & religiosity_dummy==0]
treatment_Sunni.index.nr <- support_general[treatment==2 & religiosity_dummy==0]
treatment_Muslim.index.nr <- support_general[treatment==3 & religiosity_dummy==0]
treatment_Cost.index.nr <- support_general[treatment==4 & religiosity_dummy==0]
treatment_Sunni_Cost.index.nr <- support_general[treatment==5 & religiosity_dummy==0]
treatment_Muslim_Cost.index.nr <- support_general[treatment==6 & religiosity_dummy==0]

groups <- seq(1,6,1)

groups.index.nr <- cbind( control.index.nr, treatment_Sunni.index.nr, treatment_Muslim.index.nr, treatment_Cost.index.nr, treatment_Sunni_Cost.index.nr, treatment_Muslim_Cost.index.nr)

means2 <- colMeans(groups.index.nr, na.rm=TRUE)

se.mean <- function(x){
  x.nona <- x[!is.na(x)]
  N <- length(x.nona)
  sd <- sqrt(var(x.nona))
  se <- sd/sqrt(N)
  return(se)
}


ses2 <- apply(groups.index.nr, 2, se.mean)

uis2 <- means2 + 1.96*ses2
lis2 <- means2 - 1.96*ses2


means <- c(means1,means2)
means <- means[c(1, 7, 2,8 ,3, 9, 4, 10, 5, 11, 6, 12)]

ses <- c(ses1, ses2)
ses <- ses[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

uis <- c(uis1, uis2)
uis <- uis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]

lis <- c(lis1, lis2)
lis <- lis[c(1,7,2,8,3, 9, 4, 10, 5, 11, 6, 12)]


data6.r <- as.data.frame(cbind(x1,condition, means, ses, uis, lis, group),
                         use.value.labels= T)


P6.r <- ggplot(data6.r, aes(x=x1, y=means,colour=factor(group))) +      
  geom_errorbar(aes(ymin=lis, ymax=uis,width=0.2),
                size=.5,    # Thinner lines
                width=.3,
                position=position_dodge(.9)) +
  geom_point(aes(colour=factor(group),shape=factor(group)),size=4)+
  scale_colour_discrete(name  ="",
                        breaks=c(1,2),
                        labels=c("Highly Religious", "Less Religious"))+
  scale_shape_discrete(name  ="",
                       breaks=c(1,2),
                       labels=c("Highly Religious", "Less Religious"))+
  geom_vline(xintercept=c(2.5,4.5,6.5,8.5,10.5),color="white", lty=1,size=4)+
  ylab("Mean response") +
  ggtitle("Figure 4.6. Composite Index\n") +
  coord_cartesian(ylim=c(-1, 1))+
  theme(axis.title.x = element_blank(),legend.position=c(0.9,1))+
  scale_x_continuous(breaks=c(1.5,3.5,5.5,7.5,9.5,11.5),labels=c("Control \n", "Sunni \n", "Muslim \n", "Cost \n", "Sunni\nx\ncost", "Muslim\nx\ncost"))

multiplot(P1.r, P4.r, P2.r, P5.r, P3.r, P6.r, cols =3) 





